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V  ABSTRACT 

Simultaneous  autoregressive  models  are  applied  to  balance  incomplete 
block  and  lattice  square  designs.  Relationships  between  the  resulting  maximum 
likelihood  estimators  and  standard  analyses  are  explored  and  clarified.  Also, 
it  is  shown  that  the  Papadakis  estimator,  as  applied  to  these  designs,  has  the 
characteristics  of  the  standard  analyses,  but  simply  involves  fewer  plots  in 
the  adjustment  procedure  for  treatment  averages.  A  worked  example  illustrates 
the  differences  in  treatment  estimates  in  a  lattice  square  design,  for 
standard  inter-block  analysis  and  for  simultaneous  autoregressive  assumptions. 
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SIGNIFICANCE  AND  EXPLANATION 


The  problem  of  analyzing  data  from  field  trials  when  observations  in 
neighboring  or  nearby  plots  are  correlated  is  an  important  one  which  has  been 
extensively  examined  over  the  years.  In  1937,  Papadakis  suggested  a  sensible 
but  mysterious  treatment  estimator  form  whose  role,  and  relationship  to  the 
corresponding  maximum  likelihood  estimator  has  been  debated  for  forty-odd 
years.  These  matters  were  clarified  for  certain  one-  and  two-dimensional 
designs  under  certain  simultaneous  autoregressive  models  by  Draper  and  Faraggi 
(1984).  In  the  present  paper,  simultaneous  autoregressive  models  are  applied 
to  balanced  incomplete  block  and  lattice  square  designs.  Relationships 
between  the  resulting  maximum  likelihood  estimators  and  standard  analyses  are 
explored  and  clarified.  Also,  it  is  shown  that  the  Papadakis  estimator,  as 
applied  to  these  designs,  has  the  characteristics  of  the  standard  analyses, 
but  simply  involves  fewer  plots  in  the  adjustment  procedure  for  treatment 
averages.  A  worked  example  illustrates  the  differences  in  treatment  estimates 
in  a  lattice  square  design,  for  standard  inter-block  analysis  and  for 


The  responsibility  for  the  wording  and  views  expressed  in  this  descriptive 
summary  lies  with  MRC,  and  not  with  the  authors  of  this  report. 


USE  OF  SIMULTANEOUS  AUTOREGRESSIVE  MODELS  IN  BALANCED 
INCOMPLETE  BLOCK  AND  LATTICE  SQUARE  DESIGNS 

N.  R.  Draper  and  D.  Faraggi 

1.  INTRODUCTION 

The  analysis  of  field  trials  when  observations  in  nearby  or  neighbouring 

plots  are  correlated  has  received  extensive  study  over  a  nuirber  of  years. 

An  early  reference  point  is  the  work  of  Papadakis  (1937).  The  relationship 

of  his  suggested  estimator  to  the  maximum  likelihood  estimators  in  certain 

one-  and  two-dimensional  designs  under  certain  simultaneous  autoregressive 

models  has  been  explored  by  Draper  and  Faraggi  (1984). 

It  has  been  suggested  by  Kempton  and  Howes  (1981,  p.  65  and  p.  69) 

that  the  Papadakis  model  may  be  used  quite  generally  for  field  trials  and 

specifically,  for  example,  for  lattice  squares  (LS's).  On  the  other  hand,  a 

2 

standard  analysis  for  such  designs  under  the  "usual  N(0,  Io  )"  type  of  error 
assumption  is  due  to  Yates  (1936,  1939,  1940).  Under  Yates'  analysis,  ad¬ 
justments  are  made  to  the  treatment  averages  for  the  incompleteness  of  the 
blocks. 

In  this  paper,  we  investigate  lattice  squares  under  the  assumption  that 
a  simultaneous  autoregressive  model  which  assumes  correlation  within  rows 
and  columns  is  appropriate.  This  is  done  with  different  correlations, 

Pr  for  rows  and  p£  for  columns.  We  then  obtain  maximum  likelihood  estimators 
for  treatments  and  correlations  and  discuss  the  relationships  and  differences 
between  the  three  possible  treatments  estimators  namely  maximum  likelihood 
estimators,  Yates',  and  Papadakis.  It  will  be  shown  that,  for  these  designs. 

Sponsored  by  the  United  States  Amy  under  Contract  No.  DAAG29-80-C-0041  and 
the  Wisconsin  Alumni  Research  Foundation  via  the  UW  Graduate  School. 
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the  Yates  estimator  Is  very  close  to  the  first  order  term  In  an  expansion, 
where  valid,  of  the  maximum  likelihood  estimator.  Moreover,  the  Yates' 
and  Papadakls  corrections  to  the  treatment  means  have  similar  forms.  However, 
neither  is  actually  needed  in  designs  completely  balanced  for  rows  and  columns 

because  an  explicit  solution  can  be  obtained.  An  example 
is  presented  for  illustration  in  Section  5. 

The  lattice  square  is  a  two-dimensional  design.  As  a  preliminaryt  we 
investigate  what  can  be  regarded  as  a  one-dimensional  equivalent,  namely 
the  balanced  incomplete  block  design.  Results  and  conclusions  which  closely 
parallel  those  of  the  lattice  square  are  obtained. 
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2.  BALANCED  INCOMPLETE  BLOCK  DESIGNS 


Suppose  we  wish  to  examine  t  treatments  in  b  blocks  of  size 
k  <  t.  Suppose,  further,  that  each  treatment  occurs  the  same  nunber,  m, 
of  times  and  appears  X  times  with  every  other  treatment.  Then  the  design 
is  said  to  be  a  balanced  incomplete  block  design  (BIBO).  Note  that  there 
are  n  =  bk  =  tm  plots,  and  that  X(t-l)  =  m(k-l).  For  the  standard  details 
of  analysis  reproduced  below  see,  for  example,  Cochran  and  Cox  (1957),  or 
Yates  (1940).  The  usual  model  considered  is 

y.  =  u  +  a*  +  B.  +  e.  (2.1) 

■^js  M  s  pj  js 

if  "replicates”  are  not  a  factor;  if  they  are,  some  minor  differences  occur 
in  some  formulas.  Here  y.  is  the  observation  on  treatment  s  in  block  j, 
w  is  an  overall  mean,  a*  is  the  effect  of  the  sth  treatment,  B.  is  the 

o  J 

effect  of  the  jth  block  and  e.  is  an  error  such  that  e  =  (column  vector  of 

2 

e.  in  a  defined  order)  ^  N(0,lo  ).  Two  alternative  assumptions  on  the  Bi 
js  j 

are  commonly  made. 

(a)  B.  is  a  fixed  effect.  This  leads  to  the  so-called  intra-block 

J 

analysis. 

(b)  B  =  (B, ,. . .  ,B,J '  ~  N(0,lo?)  independently  of  the  e.c.  This  leads 

I  D  P  J5 

2 

to  the  analysis  with  the  recovery  of  inter-block  information.  As  Og  °°, 
case  (b)  -*■  case  (a). 


We  set 
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Bg  a  total  of  all  response  observations  y  from 
blocks  in  which  treatment  s  appears. 

T$  *  total  of  all  observations  on  treatment  s. 

G  =  Z  T  =  grand  total  (=  I  EL  ). 
s  s 

Then,  from  Cochran  and  Cox  (1957,  pp.  445-446)  we  can  estimate  the  effect 
of  the  sth  treatment  by 

a*  =  m‘1{Ts+0[(t-k)Ts-(t-l)Bs+(k-l)G]}  (2.2) 

/v 

where  0  is  an  adjustment  factor  taking  the  value 

0  =  (t(k-l)}"1  (2.3) 

in  case  (a)  when  inter-block  information  is  not  recovered,  and  the  value 

0  =  (b-l)(Eb-Ee)/{t(k-l)(b-l)Eb+(t-k)(b-t)Ee>  (2.4) 

where  it  is  recovered.  The  quantities  Eb  and  Eg  are  the  mean  squares  due 
to  blocks  adjusted  for  treatments,  and  residual,  respectively.  Define  the 
t  by  bk  =  n  matrix  T  as  follows.  When  the  blocked  design  is  listed  in 
order,  block  by  block,  T  ^  =  1  if  the  sth  treatment  is  associated  with  the 
ith  observation,  and  T  ^  =  0  otherwise.  Let 


D 


k 


IT 


(2.5) 


Si 
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whene  T  *  (1,1,...,1)  of  dimension  1  by  k.  Then  define 

D  =  Ib  ©  Dk  (2.6) 

where  ©  denotes  the  Kronecker  product.  Then,  in  matrix  form,  with 
a*  =  (at,. . . ,a*,. .. ,a|) ' ,  (2.2)  becomes 

a*  *  m~^{TY-e(t-l)TD[I-m"^T'T]Y}  (2.7) 

where  Y  is  the  vector  of  observations,  recorded  in  the  same  pattern  as  the 
block  by  block  enumeration. 

We  now  consider,  instead  of  case  (a)  or  case  (b)  above,  an  alternative 
model  y.j  *  °s  +  »  where 

*i  *  0  j|,  "ijxj  +  Ei  <2-8> 

and  where  n^  *  1  if  y^  and  yj  occur  In  the  same  block,  and  n^j  *  0  otherwise. 
We  assume  c  *  (c^,...»cn)'  -  N(0,lc  )  independently  of  the  x's.  Following 
Besag's  (1974,  Eq.  (4.13)  with  y  =  0)  simultaneous  autoregressive  model, 
we  obtain  the  likelihood 


(2TO2)"n/2|B|exp{-(2a2)"1X,B,BX},  (2.9) 

A 

where  B  »  I  -  pD  and  X  *  Y  -  T'a.  The  maximum  likelihood  estimator  a 
is  of  form  (Ripley,  1981,  (5.34)) 
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(TV"1T')~1TV”1Y1  (2.10) 

where  TV_1T'  =  T(I-pD)2T'  =  m{I-m”V(2pD-p2D2)T'}  using  the  fact  that  TT'  =  ml. 
The  matrix  (TV"^T')~^  can  be  expanded  as  a  convergent  series  if  and  only 
if  the  eigenvalues  of  Qg  =  m  T(2pD-p  D  )T'  are  all  less  than  one  in  absolute 
value.  This  is  true  if 

r  min{(l+21/2),t1/2-l),  k  =  2, 

(l-21/2)/(k-l)  <  p  <  )  min{^(l+21/2),  J(t-1)},  k  =  3,  (2.11) 

L  (l+21/2)/(k-l),  k>4. 

If  we  retain  only  terms  of  orders  zero  and  one  in  Qg,  we  obtain 

*  m"1[TY-{2p-pZ(k-2)}TD(I-m“1T'T)Y]  (2.12) 

where  we  have  used  the  facts  that  D2  =  (k-l)In  +  (k-2)D,  and  T(I-m”V'T)  =  0. 
Comparing  (2.12)  with  (2.7)  we  notice  the  estimators  are  identical  except 
that  e(t-l)  in  (2.7)  replaces  (2p-p2( k-2) >  In  (2.12).  Because  of  the  validity 

A  A 

of  the  expansion,  will  converge  to  the  maximum  likelihood  estimator  a, 

A 

for  fixed  p,  if  the  former  is  applied  iteratively  in  the  following  sense: 

aj+1  *  m”^[TY-{2p-p2(k-2)}TD(Y-T'oij)]  (2.13) 

*  -1 

where  aQ  =  m  TY.  In  fact,  iteration  is  not  needed  here  because  the  balanced 
incomplete  block  design  produces  a  TV”^T'  of  patterned  form  which  can  be 
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explicitly  inverted  (Rao,  1973,  p.  67)  to  give  o  in  terms  of  p  as 

o  =  a( I+bTDT* )T( I-Jd)Y  (2.14) 

where  a  =  {m-(t-2)<J>}/[{m-(t-lHHm+(|>}],  b  =  4>/(nr:'( t-2)<p}  and 

A  /V  AO  AO  A 

$  a  {2p-(k-2)p  }/{l+(k-l)p£).  The  maximum  likelihood  p  is  obtained  as 
follows.  Setting  the  first  derivative  of  the  log  likelihood  with  respect 
to  p  equal  to  zero  provides 


(Y-T'ct) 1  (D-pD;  (V-T'a)  .  (k-l)p  (2  15 

(Y-T  a)(I-2pD+p^D  ;(Y-T'a)  {l+(n-l)p}{l+^} 

A  A 

The  a  may  be  substituted  from  (2.14)  to  provide  an  equation  in  p  which 
can  be  solved,  for  example,  via  the  Newton- Raphson  method.  From  this,  a 
is  then  obtained.  In  practice,  it  is  easier  to  replace  a  in  the  likelihood 

A 

by  a  from  (2.14)  and  then  to  maximize  the  likelihood  numerically  by  a 

A 

search  over  p. 

For  least  squares  estimation,  the  term  on  the  right  of  (2.15)  is 
replaced  by  zero,  so  that 


PLS  *  (Y-TaLS)'D{Y-T^LS)/{(Y-T^LS),D2(Y-TaLS)>,  (2.16) 


A  A  A  A 

where  has  the  same  form  as  a  but  with  p  replaced  by  p^.  Again 
numerical  maximization  is  simplest. 


I 


-8- 

3.  LATTICE  SQUARES 

Suppose  we  wish  to  examine  t  treatments  in  a  set  of  m  k  by  k 
lattice  squares.  Typically  m  =  6(k+l)  where  5  =  j-,  1,  2,  ...  .  The  different 
values  of  6  correspond  to  various  association  relationships  among  the 
treatments.  When  5  =  any  treatment  s  will  appear  with  any  other  treatment 
w  once  in  either  a  row  or  a  column  but  not  both.  When  6  =  1,  s  and  w 
appear  together  both  in  a  row  and  a  column,  once  each,  and  so  on.  The  usual 
model  considered  is 

yi  =  u  +  as  +  *q  +  &r  +  Yc  +  ei  (3-1) 

where  the  set  of  subscripts  (q,r,c)  representing  the  position  in  the  rth 
row  and  the  cth  column  within  the  qth  replicate  have  been  replaced  by  a  single 
subscript  i  =  l,2,...,n  which  denotes  the  plot's  position  in  a  row  by  row 
enumeration  of  the  plots,  replicate  by  replicate.  Two  alternative  assumptions 
on  the  B's  and  the  y's  are  commonly  made 

(a)  Br»  Yc  are  fixed  effects.  This  leads  to  an  intra-block  analysis. 

2  2 

(b)  8_  ..  N(0,cC)  and  v  ~  N(0,a  )  and  all  B's  and  y's  are  independent. 

r  p  c  Y 

2  2 

Typically  it  is  assumed  that  =  ag,  although  this  is  not  necessary.  As 
2  2 

oa  00  and  a  ■*  ®,  case  (b)  -*■  case  (a). 

P  Y 

He  set,  for  6  >  1 , 
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T$  =  total  of  all  response  observations  y  from  all  plots  i 
which  receive  treatment  s, 

Rs  =  total  of  all  y's  in  rows  in  which  treatment  s  appears, 

C$  =  total  of  all  y's  in  columns  in  which  treatment  s  appears 

G  =  grand  total  of  all  y's, 

Ls  *  (m-l)Ts  -  mR$  +  G, 

M$  =  {m-l)Ts  -  mC$  +  G. 


> 


Then,  from  Cochran  and  Cox  (1957,  p.  493),  we  can  estimate  the  effect  of 
the  sth  treatment  by 


a*  =  m'1{Ts+XrLs+X(.Ms},  (3.2) 

where  Xr  and  Xc  are  coefficients  for  the  adjustments  needed  for  the  in¬ 
completeness  of  the  rows  and  columns.  The  values  of  Xr  and  Xc  depend  on 
the  replication  selected  and  which  case,  (a)  or  (b)  is  assumed.  Examples 
are  in  Appendix  A. 

For  6  *  £■,  the  formula  for  a*  is  of  the  same  form  as  (3.2)  but  with 
redefined  symbols.  We  do  not  discuss  this  because  our  matrix  formulation 
does  not  contain  this  ambiguity. 

Define  the  t  by  n  matrix  T  as  follows.  When  the  blocked  design 
is  written  out  row  by  row  within  replicate  and  then  replicate  by  replicate, 
T$^  *  1  if  the  sth  treatment  is  associated  with  the  1 th  observation,  and 
Tsi  =  0  otherw'ise-  Let  =  11*  -  1^,  where  T  =  ( 1 .1 ....  ,1 )  of  dimension 
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1  by  k,  and  define 


2r  ‘  ©  ]lc  ©  V 

zc  *  >m  ©  °k  ©  V  <3-3> 


where  (x)  denotes  the  Kronecker  product.  Then,  In  matrix  form,  with 
oL  *  (cij ,<*2 ♦  •  •  •  »at) '  *  tf,e  *■  <fenot1ng  lattice,  (3.2)  becomes 


aL  -  m'1{TY-mT(XrZr+XcZc)(I-m^1T,T)Y} 


(3.4) 


where  Y  Is  the  vector  of  observations,  recorded  In  the  same  pattern  as 
described  earlier. 

We  now  consider.  Instead  of  case  (a)  or  case  (b),  an  alternative 
model  y^  3  “s  +  x^»  1  3  1,...,n,  s  *  l,...,t,  where 


o-  £  n(^x.  +  P -  £  ni^x. 
r  j+1  j  c  J+1  j 


(3.5) 


and  where  n!^  *  1  If  y.,  and  y1  occur  in  the  same  row  and  replicate  and  is  zero 
J  Ic)  J 

otherwise,  while  -1  If  y^  and  y^  occur  In  the  same  column  and  replicate  and 

2 

Is  zero  otherwise.  We  assume  e  ~  (e^,...,en)'  ~  N(0,lo  )  independently  of  the  x' 
Following  Besag's  (1974,  Eq.  (4.13)  with  y*0)  simultaneous  autogressive  model,  we 
obtain  the  likelihood  (27to2)'n/2|B|exp{-(2o2)"1X 'B'BXl  where  B=I-prZr-pcZc, 


1 


-n- 


X  =  Y  -  T’o. 
where  TV’ W* 


*  _i  _i 

The  maximum  likelihood  estimator  a  is  of  form  (TV  T ' )  TV  Y, 
=  T(I-prZr-pcZc)^T'  =  m(I-QL)  and  where 

-  m-1T(2^rZr+2JcZc-^Z2.;2z2_2— cZrZc)T. .  (3.6) 


The  matrix  (TV’V)"^  can  be  expanded  as  a  convergent  series  if  and  only  if 
the  eigenvalues  of  QL  are  all  less  then  one  in  absolute  value.  This  is  true 
if 


-1  <  q1  +  (k2-l )q2  <  1, 

-1  <  ^  -  q2  <  1 . 


where 


q1  -  -(m-2)(p*+p*), 

q2  -  {2pr+2pc-(m-3)(p^+p^)  -  2(k-l)prpc}/(k+l).  (3.7) 


For  example,  if  m  *  6,  k  *  5,  the  conditions  imply  a  feasible  region  shown 
shaded  In  Figure  1. 

If  we  retain  only  terms  of  orders  zero  and  one  in  Q^,  we  obtain 


-  m’1[TY-T(JrZr4cZc)(I-  m^T'TjY], 


(3.8) 


Figure  1.  The  shaded  region  shows  the  values  of  (pr,pc)  for  which  the  ex¬ 
pansion  of  the  maximum  likelihood  estimator  Involving  QL  is  valid 
for  the  case  m  *  6,  k  *  5  described  in  the  text.  For  the  case 
Pr  =  Pc  =  P  thl's  region  reduces  to  the  intersection  of  the  shaded 
region  and  the  45°  line  shown,  with  a  multiplicative  scaling 
factor  of  2”^  applied. 
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A  A  A0  /\  A  A I 

where  4>u  *  2py  -  (m-3)p‘  +  2p rp£  and  u  *  r,  c  In  turn.  In  deriving 
we  have  used  the  following  facts  (u=r  or  c): 


ZZ  *  (m-2) I  +  (m-3)Zu, 

TZrZc  -  1t1iT-T-TZr-TZc 

TU-nfVT)  =  0, 

I-m"”1  T *T)  *  0. 


(3.9) 


Comparing  (3.8)  with  (3.4),  we  see  that  the  estimators  are  identical  except 

A 

that  mXu  in  (3.4)  replaces  4>u  in  (3.8).  Because  of  the  validity  of  the 

A..  . 

expansion,  a^-  will  converge  to  the  maximum  likelihood  estimator  a  ,  for 

A  A  Aj 

fixed  pr  and  pc,  if  Is  applied  iteratively  in  the  following  sense 

?j+1  *  m’1{TY-T(JrZr4cZc)(Y-T1^)},  (3.10) 

-I 

where  aQ  »  m  TY.  In  fact,  iteration  is  not  needed  here  for  6  1,  because 

the  lattice  design  produces  a  TV'H'  of  patterned  form  which  can  be  explicitly 

Al  A  A 

inverted  (Rao,  1973,  p.  67)  to  give  a  in  terms  of  p  and  p  as 

•  v 


((a-b)I+bll '}T(I-Q)Y, 


(3.11) 


where 


a  =  {(m-l)q1-(k2-2)nxi2}/C{(m-1)q1-(k2-1)mq2H(m-l)q1+mq2}], 
b  =  mq2/[{(m-l)q1-(k2-l)mq2}{(m-l)q1+mq2)].  (3.12) 


The  maximum  likelihood  pp  and  pc  are  obtained  as  follow.  Setting  the  first 
derivatives  of  the  log  likelihoods  with  respect  to  p^  and  p2  equal  to  zero 
provides  equations: 


nX'ZuBX 

A  A^A 

X'B  a 


fe*  |8‘]  -  * 

L  u  Jp  =p  ,p  =p 

- r. . r— =  o, 

I B  | 


r,c 


(3.13) 


where  X  and  B  are  the  values  of  X  and  B  with  a  ,  p^  and  pc  replacing 

a,  pr  and  pc<  These  two  simultaneous  equations  involve,  in  their  second  terms, 

the  (1,1+1)  and  the  (1,1+k)  elements,  the  designations  being  reduced  modulo 

A.1 

n  where  they  exceed  n,  respectively  of  B  and  must  be  solved,  a  tedious 
calculation.  It  is  simpler  in  practice  to  substitute  a  selected  grid  of 
(pr.Pc)  values  into  oL  and  hence  evaluate  the  likelihood  over  that  grid,  so 
picking  out  the  maximum  value  of  the  likelihood. 

For  least  squares  estimation,  the  term  on  the  right  of  (3.13)  is 


replaced  by  zero  so  that 
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“  A  " 

pr,LS 

s 

A 

pc,LS 

rr 


rc 


cc 


■  -1 

X'ZrX 

X'Z  X 

c 

(3.14) 


where  buy  =  X'Z^ZyX.  The  least  squares  estimate  of  a,  a^s  is  of  the  same 

A I  A  A 

form  as  a  but  has  p  replaced  by  p  .  everywhere. 

U  U  )L j 

Special  Case  pp  *  pc  *  p. 

In  many  practical  situations,  it  is  reasonable  to  assume  that 
Pr  =  Pc  =  P*  For  comments  on  this  choice,  see  Kempton  and  Howes  (1981,  p.  63). 
All  formulas  reduce  to  their  appropriate  forms  via  simple  substitution  of 
Pr  *  pc  *  p-  All  that  is  needed  is  to  write  Z  =  Zp  +  Zc  whenever  it  occurs 
and  to  drop  subscripts  from  A,  <p,  and  p.  The  equation  parallel  to 
Eqs.  (3.13)  is  the  result  of  adding  these  two  equations  together. 


1 

i 
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4.  LATTICE  SQUARE  DESIGN  EXAMPLE 


We  shall  re-examine  data  given  by  Cochran  and  Cox  (1957,  pp.  490-493) 
comprising  a  lattice  square  design  with  k  =  4,  6  -  1 ,  m  =  5,  n  =  80.  The 
appropriate  formulas  for  Yates'  estimators  Ar  and  Ac  are  given  in  Appendix  A 
and  here  result  in  the  numerical  values  A  =  0.04787  and  A  =  0.03037  (p.  493). 

r  v 

/V  A 

The  maximum  likelihood  estimates  for  p„  and  p„  are  p„  =  0.108  and  p  =  -0.030 

r  c  r  c 

and  these  lead  to  estimates  for  the  16  treatment  parameters  given  in  the 
second  column  of  Table  1.  These  may  be  compared  with  the  "adjusted  means" 
of  Cochran  and  Cox  (1957,  p.  491)  which  are  the  corresponding  Yates'  estimates, 
shown  as  the  third  column  of  Table  1.  We  see  that,  although  two  entirely 
different  models  are  used,  individual  estimates  are  broadly  speaking,  compatible. 
In  view  of  the  fact  that  only  one  at  most  of  these  models  can  be  correct  for 
the  data  set,  we  do  not  compare  them  via  a  model -dependent  criterion  such  as 
weighted  or  unweighted  sum  of  squares  of  residuals. 


Table  1.  Comparison  of  treatment  effects  for  a  4  by  4  lattice  square  design 
with  five  replications.  Second  column:  Maximum  likelihood  estimates  under 
the  assumption  of  a  simultaneous  autoregressive  model.  Third  column:  Yates' 
estimates  from  Cochran  and  Cox  (1957,  p.  491). 


Treatment 

Mle 

Yates' 

1 

4.98 

6.45 

2 

12.71 

13.68 

3 

8.90 

8.73 

4 

11.41 

11.36 

5 

9.89 

9.44 

6 

5.96 

7.58 

7 

7.52 

7.37 

8 

9.79 

9.32 

9 

10.91 

10.01 

10 

15.77 

14.91 

11 

18.12 

17.59 

12 

12.91 

12.70 

13 

12.07 

10.69 

14 

13.40 

14.27 

15 

10.07 

9.28 

16 

10.05 

11.09 

Sum 

174.46* 

174.47* 

♦These  figures  are  subject  to  small  rounding  errors; 
both  should  be  174.48. 
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5.  PAPADAKIS  AND  MAXIMUM  LIKELIHOOD 

We  reconsider  the  assumption  made  in  connection  with  Eq.  (3.1). 

Suppose  that  now  $r  is  the  combined  neighbouring  effect  from  the  "plot  to 
the  left  In  a  row",  and  the  "plot  to  the  right  In  a  row",  these  left  and 
right  adjacencies  being  determined  In  torus,  wrap-around  manner  (Martin, 
1982).  A  similar  "up  and  down  in  a  coluim"  definition  Is  applied  now  to 
Yc.  Let  Nr  and  Nc  be  the  n  by  n  neighbour-specification  matrices  for  rows 
and  columns  whose  ith  row  contains  1  in  positions  j  for  which  plot  j  is 
row-  and  column-adjacent  respectively,  to  plot  1.  and  zero  otherwise.  Then 
the  maximum  likelihood  estimator  of  a*  is 

a*  =  m”1{TY-mT(^r/ir^cNc)(I-m‘1T,T)Y)  (5.1) 

where  ij>r  and  4<c  are  coefficients  for  the  adjustments  needed  for  the  row- 
and  column -adjacency  effects.  Parallel  to  the  work  in  Section  3,  there  will 
be  two  forms  for  and  ^  according  to  whether  an  intra-neighbour  analysis 

•  V* 

(BrtYc  fixed  effects)  or  an  Inter-neighbour  analysis  (Br,Yc  random  effects) 
is  specified. 

If  we  compare  (5.1)  with  Eq.  (3)  of  Draper  and  Faraggi  (1984),  we 

A. 

see  that  they  would  be  equivalent  If  n*J>r  *  tnpc  ■  $  and  if  N  *  Nr+Nc*  *n 
fact,  (5.1)  Is  the  Papadakis  estimator  for  the  case  f  \|»c,  up  to  definition 

of  and  i|;  .  The  Papadakis  formulation  would  Involve  regression  estimators 
of  and  4>c  Instead.  We  thus  see  that  the  Papadakis  estimator  can  be 
regarded,  up  to  definition  of  il>r  and  i|/  ,  as  a  maximum  likelihood  estimator 
In  the  particular  model  formulation  described  above. 
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6.  REMARKS 


A  previous  paper  (Draper  and  Faraggi,  1984)  discussed  field  trials  in 
one  and  two  dimensions  for  designs  of  Type  11(a),  Type  III,  and  Type  III 
cyclically  row  permuted  (Williams,  1952)  for  one  and  two  dimensional  layouts. 
In  the  two  dimensional  layout,  no  directional  differences  were  assumed  and  a 
single  p  value  was  used  for  both  rows  and  columns.  The  extension  to  pp  and 
Pc  for  that  case  Is  easily  effected  using  formulas  of  the  type  exhibited  in 
this  paper.  The  maximum  likelihood  estimator  has  the  form  of  (2.10),  with 

V*1  -  (I-PrNr-PcNc)2  (6-D 


where  Nr  Is  an  c  by  n  row-nelghbour-speclflcatlon  matrix  defined  by 


Nr  *  I«(k+l)©Ik©N’ 


(6.2) 


where  Nc  Is  an  n  by  n  column-neighbour-specification  matrix  defined  by 


Nc  *  'S(W)©N®Ik- 


(6.3) 


and  where  N  Is  the  cyclic  k  by  k  matrix 


0  1  0  0  ...  0  1 

1  0  1  0  ...  0  0 

•  •  •  • 

1  0  0  0  1  0 


(6.4) 


We  can  re-write  (2.10)  as 


m‘1(I-fn-1TQT,)‘1T(I-Q)Y 


(6.5) 


where 

Q  ■  2jrNr  *  2?CNC  -  SjjlJ  -  7/c  -  !PrPcHr1c-  (6.6) 

Here, 

N2  »  21  +  n[.2) ,  N2  =  21  +  N<2) ,  NrNc  “  N^1 } , 

where  N^2^ ,  N^2^ ,  and  denote  n  by  n  neighbour-specification  matrices 

i  c  re 

which  contain  I's  and  0's  such  that 


"r2)  ’  lS(W)©Ik©N'2) 

(6.7) 

"£2)  *  ■S(k+1)©N‘2)©Ik 

(6.8) 

"t'c  1  ■  'S(W)GN©N 

(6.9) 

where  Is  the  cyclic  k  *  k  matrix 

/ 


N'  '  = 


0  0  10 


•  •  • 


0  1  0 


0001  ...001 

1  0  0  0  1  ...  0  0  0 
•  •  •  • 

0100  ...100 


(6.10) 


and  N  Is  defined  In  (6.4).  For  example,  the  ith  row  of  has  Vs  in 
those  columns  j  for  which  cell  j  is  located  next  but  one  to  cell  i 
In  a  row  of  the  design,  and  has  zeros  otherwise.  Table  2  illustrates  the 
general  situation  and  has  those  (row)  cells  marked  with  symbol  N^K 
Similar  remarks  apply  for  and  for  columns  and  diagonals 

respectively,  as  shown  in  Table  2.  When  cell  i  is  near  or  on  an  edge,  the 
appropriate  torus-generated  cells  must  be  used. 


N«> 

rc 

Nc 

„<*) 

Nr 

Cell 

1 

"r 

n(2) 

N<"> 

rc 

Nc 

rc 

Table  2.  Cell  patterns  which  lead  to  the  n  by  n  neighbour-specification 
matrices  Nr,  Nc,  nJ.2),  N<2), 
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Appendix  A.  Values  of  Xr,  Xc  for  i  =  1,  2  for 

Lattice  Square  Designs. 


6 

With  recovery  of 

Inter-block  Information 

Without  recovery  of 
Inter-block  Information 

xr 

Xc 

xr  *  xc 

1 

2(Er-Ee) 

2(EC-Ee) 

2 

2 

k(k+l)Er 

k(k+l)Ec 

k(k+l ) 

1 

(Er“Ee)(kEc-Ee) 

<Ec-EeHkEr-Ee> 

1 

(k-l)(k2ErEc-E^) 

(k-l)(k2ErEc-E2) 

TFTTk 

2 

W  -W 
e  r 

VMc 

1 

k(wr+Wc+(k-i)we} 

k{Wr+Wc+(k-tjweJ 

TFTJk 

In  the  above  table, 

Wfi  -  Eg1,  Wr  «  (2k-l)/(2kEr-Ee),  Wc  »  (2k-l)/(2kEc-Ee), 

Er  ■  mean  square  due  to  rows  adjusted  for  treatments, 

Ec  *  mean  square  due  to  columns  adjusted  for  treatments  and  rows, 
Ee  »  residual  mean  square. 
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